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Entanglement between large numbers of quantum modes is the quintessential resource for future 
technologies such as the quantum internet. Conventionally the generation of multimode entangle- 
ment in optics requires complex layouts of beam-splitters and phase shifters in order to transform 
the input modes in to entangled modes. These networks need substantial modification for every 
new set of entangled modes to be generated. Here we report on the highly versatile and efficient 
generation of various multimode entangled states with the ability to switch between different linear 
optics networks in real time. By defining our modes to be combinations of different spatial regions 
of one beam, we may use just one pair of multi-pixel detectors each with M photodiodes in order 
to measure N entangled modes, with a maximum number of N=M modes. We program virtual 
networks that are fully equivalent to the physical linear optics networks they are emulating. We 
present results for N=2 up to N=8 entangled modes here, including N=2,3,4 cluster states. Our 
approach introduces flexibility and scalability to multimode entanglement, two important attributes 
that are highly sought after in state of the art devices. 
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Multi-partite entanglement is not only of fundamental 
scientific interest, it is also the key ingredient for quan- 
tum information technologies [HI]. In optics, several 
impressive demonstrations of multi-partite entanglement 
have been shown recently including an 8-photon cluster 
state jS] and a 9-mode state used for error correction 
[§]. However, these schemes tend to employ one detec- 
tion system per entangled mode/qubit, which introduces 
a lack of flexibility and is detrimental to its scalability. 
These optical setups are built to produce one set of out- 
puts or to perform one given protocol; in order to change 
the output the optical hardware itself must be modified. 
We report here on a system with the ability to switch 
in real time between desired output states using just one 
detection scheme. 

Currently the well-established recipe for generating en- 
tanglement using continuous wave laser beams is to mix 
squeezed modes of light together at beam-splitters. It 
is possible to create A-mode entanglement given a net- 
work of N-l beam-splitters with N input modes, even 
with less than N squeezed modes [7|. In our scheme we 
co-propagate all possible spatial modes of light within one 
beam. Entanglement between co-propagating modes in 
one beam has been previously demonstrated with spatial 
modes 19] , and also in the frequency domain [10] . In the 
current work we radically expand the idea of one-beam 
entanglement by introducing the notion of emulating lin- 
ear optics networks, by programming virtual networks 
that mix together different spatial regions of the light 
beam. These software based networks calculate the pre- 
cise weighted combinations of the spatial regions required 
to emulate the physical networks. This is possible be- 
cause the linear optical components in a typical network 
simply perform reversible operations, and can be repre- 
sented by unitary matrices. It is worth stating explicitly 
that the entangled spatial modes that we produce are 
event-ready and unconditional before the detection pro- 
cess. The real-time virtual networks allow us to match 
the detection basis to the desired spatial mode basis con- 
tained within the beam, analogous to shaping a reference 
local oscillator beam. 

Figure [T] shows two such networks that produce a 
2-mode entangled state and an 8-mode entangled state 
respectively. We also program virtual networks for 3, 
4, 5, 6, and 7- mode entangled states, the results of 
which are shown in Table 1 and Figure [4] As a further 
demonstration of the versatility of our setup we produce 
linear 2, 3, and 4 mode cluster states, which are highly 
entangled graph states garnering attention for their 
potential in quantum computing [2 . Table 2 summarises 
the cluster state measurements while Figure [5] shows 
that cluster states demand more stringent squeezing 
requirements than non-cluster inseparable states. 

Results 

Measuring spatial modes. By employing custom 



made multi-photodiode-homodyne-detectors (MPHD) 
that each contain an array of 8 photodiodes (see Fig- 
ure^ we detect the light in 8 spatial regions and assign 
individual electronic gains to each spatial region. The 
linear combination of the 8 gain-adjusted photocurrents 
constitutes the measurement of one mode. 

More generally, we can express the measurement pro- 
cess of a complete set of spatial modes in one beam by 
the following: 

a = Ui (1) 
= U nct U™i (2) 

where a = (a%,--- ,o,n) t is the set of N measured 
modes projected by the Nx8 unitary matrix U acting 
on the 8 homodyne-subtracted photocurrent operators 
i = (ii, ■ ■ ■ ,is) T . C/£ is an N x 8 matrix made up of 
the top N rows of U- m) the orthogonal 8x8 unitary ma- 
trix that recovers the important set of 8 unmixed spatial 
modes that span the input basis (see Methods). Input 
modes are then mixed via U net , which emulate linear op- 
tics networks, given by the N x N matrix: 

/ v\ v\ ■ ■ ■ vjf \ 
vf vi, ■■■ v z N 
C/„ot= . . . " , (3) 

v< ^ ■■■ <) 

where «™ 6 R. 

This allows us to uniquely define a mode d n by 
the 8 real numbers in the nth row of U, which we 
will label as the mode's gain vector G„, such that 
«n = G n i. Therefore each spatial mode we measure, 
whether belonging to the input basis or an entangled 
mode basis, is defined by a unique pattern within 
the light beam. These spatial mode patterns, repre- 
sented by Gaussian profiles modulated by respective 
electronic gains G„, are shown visually in Figure [3] 
while the detection stage of Figure [2] shows how we 
implement this experimentally. The spatial modes are 
orthogonal to each other, spanning a basis so that the in- 
dependent measurement of each mode is possible [TT1I12) . 

Input basis. We create two amplitude squeezed 
modes via optical parametric amplification (OPA). The 
first mode is converted to a flip mode (FM) by phase de- 
laying half its beam by half a wavelength, n (see inset of 
Figure [2]). The FM is overlapped in quadrature with the 
Gaussian mode (GM) output of the second OPA upon 
reflection of its output coupler [13] • These two squeezed 
modes are the first two modes of what we refer to as the 
input basis; a x and a 2 . Six co-propagating vacua modes 
are measured by calculating G n vectors that are orthog- 
onal to both a\ and 0,2- These vacua modes (labelled 



3 



as-.a^) complete the input mode basis (see middle row of 
Figure [3]) . Measuring these modes amounts to matching 
the detection basis by following equation (1) and setting 
t^nct = I ( sec Methods). 

Each spatial mode is characterised by the continuous- 
variable (CV) quadrature operators x and p of the electric 
field operator. The x and p variance measurements of 
the eight modes in the input basis are shown in Fig- 
ure [§i,b. Here, ([Ax GM ] 2 ) = ([Axi] 2 ) = -4.3 ± 0.05dB 
and ([Ap FM ] 2 ) = ([Ap 2 ] 2 ) = -3.7 ± 0.05dB below the 
standard quantum noise, and the variances of the vacua 
are verified to equal quantum noise. 

Entangled mode bases. Programming a virtual net- 
work amounts to calculating the precise expression for 
U ne t ■ The unitaries we have access to in programming the 
virtual networks are beam-splitters and it phase shifts. 
U ne t is the concatenation of all of these unitaries that 
make up a linear optics network. The tt phase shift is 
equivalent to multiplying a by — 1. Note that arbitrary 
phase shifts are forbidden as each measurement natu- 
rally corresponds to detection at a fixed phase defined 
by a shared reference beam, the local oscillator. Impor- 
tantly, optimal virtual networks are calculated allowing 
for optimisation of beam-splitters due to asymmetries in 
the squeezing levels of input modes. 

The most intuitive virtual network we create is the 2- 
mode EPR state [T3] shown in Figure [T^i with 2 squeezed 
inputs. Here we engineer spatial mode patterns which 
have no spatial overlap; the left half of the beam is en- 
tangled with the right half (see the top left of Figure [3]). 
Entangled modes belonging to other bases share spatial 
overlap but are nevertheless spatially orthogonal. 

Spatial modes measured in an entangled mode basis 
are given a superscript N to distinguish them from 
modes in the input basis; a\ and a\ represent the two 
modes spanning the N=2-mode EPR basis. See the 
Methods section for details on how we create virtual 
networks for each of the N=2, 3, 4, 5, 6, 7, and 8-mode 
bases. In general we construct networks pertaining to 
N modes by concatenating N-l virtual beam-splitters 
with vacua on unused input ports. This is a highly 
efficient approach to creating multimode entanglement 
as the arduous tasks of mode matching and alignment 
are replaced with the ease of programming. 

Cluster states. Attracting attention for their po- 
tential in one-way quantum computing schemes, cluster 
states are a type of highly entangled Gaussian graph 
state [HI [16]. They satisfy the quadrature relation 
(pa —J2beN — ► 0. As infinite squeezing would require 
infinite energy and are thus unrealisable, one is limited 
to the production of approximate cluster states in the 
laboratory, and there have been demonstrations of up to 
four- mode CV cluster states thus far [T7J [TH] . 

In order to measure cluster states in one beam we 



must be able to access the correct quadratures of each 
entangled mode. Here we measure 2, 3, and 4 mode 
linear cluster states, however measuring arbitrary cluster 
shapes would require modifying the optical setup (see 
Methods). 

Discussion 

In order to verify entanglement between measured 
modes we use the well-established van Loock-Furusawa 
inseparability criteria |19j . For an A-mode entangled 
state, it is sufficient to satisfy A — 1 inseparability in- 
equalities: 

(I) 

([A(ii - x 2 )] 2 ) 

+([A(pi +p 2 + 93P3 + ■■■ + SrPn)] 2 ) < 1, 

( 4 ) 

(N-l) 

([A(x N -i -^n)] 2 ) 

+ ([A(gxpi + ... + .9N-2PN-2 +PN-1 +Pn)} 2 ) < 1. 

with free parameters <?i, to be optimised for maximum 
inseparability. We have omitted the superscript N here 
for clarity, as the above holds for any mode basis. The 
subscripts n of x and p here indicate the nth mode in 
the N-mode basis. Table 1 summarises the measured de- 
grees of inseparability for all N — 1 inequalities in each 
A-mode basis, given in roman numerals. Table 2 sum- 
marises the more stringent inseparability required for 
unweighted cluster states (all homodyne gains {gi} are 
set to 1). The relevant inseparability inequalities satisfy 
the cluster state quadrature relationship written in the 
form of equation (4) , and are written out explicitly in the 
Methods for reference. 

The terms in equation (4) measure the degree of corre- 
lations between any two modes in a given basis. For 
the modes to be inseparable each of these correlation 
variances (correlations in the x quadrature and anti- 
correlations in the p quadrature) must be in the quantum 
regime, that is below the normalised quantum noise of 
two units of vacua. Figure [I];, d shows this to be the case 
in our experimental measurements. Although we are lim- 
ited here to 8 modes due to our detection scheme, this 
scheme is scalable to higher numbers of mode entangle- 
ment even without increasing the number of squeezing 
resources, as shown in the simulation traces of Figure [5] 
As we increase the simulated number of modes in the 
basis up to 30, the degree of inseparability approaches 
the classical bound of 1 due to the vacuum noise penalty 
for each additional unsqueezed mode input. Entangle- 
ment is shown to hold here however, even with current 
squeezing levels. Importantly, there is no loss incurred 
during the transformation of the squeezed input modes 
into a set of entangled modes, as can be seen by the 
agreement of the theoretical predictions and the experi- 
mental values of Figure [5^,. This equates to perfect mode 
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matching at every virtual beam-splitter. Figure [5]d ex- 
plores how inseparability scales with different squeezing 
levels. Measuring a larger number of inseparable modes 
experimentally requires only an increase in the number 
of photodiodes in the MPHD, and importantly no modi- 
fication of the optical setup. Note that this is not true for 
cluster states, and the number of squeezed inputs must 
be increased accordingly. 

For the special case of N = 2, optimal EPR entan- 
glement [20] is measured to be 0.58 ± 0.01. Optimis- 
ing for the beam-splitter reflectivity [2H [22] we find that 
due to the slight asymmetry between input squeezing lev- 
els, the optimal beam-splitter ratio here is not 50%, but 
rather 48.8%, leading to a very slight improvement over 
the symmetric network. Each unique beam-splitter re- 
flectivity changes the mapping of U nc t such that formally 
the beam of light contains an infinite number of mode 
bases. The versatility of our scheme comes from being 
able to match the detection basis to a network that has 
been optimised for an arbitrary set of inputs. 

The entanglement demonstrated in the current work 
allows for such protocols as quantum teleportation 
[251 [21] • To perform complex protocols such as 
one-way measurement based quantum computations 
[2] modifications are needed. We have shown with the 
current setup it is possible to create cluster states, the 
resource state for one-way quantum computations. To 
perform computations on the cluster states however we 
need access to arbitrary homodyne angles of each mode, 
and the ability to perform feed-forward to any desired 
mode [5S]. Both are feasible with existing technologies 
[27] [28] as discussed in the supplementary material. 



Emulating linear optics networks by mixing coprop- 
agating spatial modes is a highly efficient method 
for generating multimode entanglement. Otherwise 
arduous and potentially lossy tasks such as mode 
matching during the construction of a linear optics 
network are performed effortlessly and losslessly via 
software controlled combinations of the spatial modes. 
We have shown that although correlations weaken if 
more squeezing resources are not added, (non-cluster) 
entangled modes scale here as the number of orthogonal 
modes measurable within the beam. The maximum 
number of measurable modes corresponds directly to the 
number of photodiodes in each pair of the multi-pixel 
detectors. We have demonstrated this by measuring 
N=2, 3, 4, 5, 6, 7, and 8-mode entanglement within one 
beam, including up to 4 mode cluster states, switching 
between them in real time. The ability to perform a wide 
range of protocols and optimise networks for asymmetry 
using just one optical setup offers versatility to future 
networks that will utilise entanglement as a resource. 

Methods 

Experimental setup. We use a dual-wavelength continuous- 



TABLE I: Inseparability of entangled modes based on 
the van Loock-Furusawa criteria. Each row shows that 
for a basis of N quantum modes, the N-l values obtained 
from quadrature variances are well below 1. This verifies 
entanglement of the N modes. 

N I II III IV V VI VII Avg. 
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0.39 












0.39 
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0.56 


0.56 










0.56 
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0.64 


0.63 


0.64 








0.64 
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0.69 


0.69 


0.70 


0.70 






0.69 
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0.73 


0.73 


0.75 


0.74 


0.74 




0.74 
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0.77 


0.78 


0.77 


0.76 


0.77 


0.77 


0.77 


8 


0.79 


0.79 


0.78 


0.81 


0.79 


0.80 0.79 


0.79 



^Uncertainty is ±0.01 in all cases. 

TABLE II: Inseparability of cluster states. 
N I II III IV 



Avg. 
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0.39 








0.39 
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0.49 


0.70 






0.59 
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0.79 


0.67 


0.84 




0.76 


5 


0.79 


0.67 


1.10 


1.18 


0.93 



''Uncertainty is ±0.01 in all cases. 



wave Nd:YAG laser at 1064 nm and 532 nm. The optical 
parametric amplifiers (OPA) each contain a periodically poled 
KTP crystal in a bow-tie cavity. The squeezed beams are almost 
identical in purity, with squeezing levels of approximately -6 dB 
and anti-squeezing of 8.5 dB. The beam containing the 8 spatially 
orthogonal modes (see main text) is made highly elliptical in 
order to be measured by the MPHD, which has a linear array of 8 
photodiodes. The photodiode array used is a Hamamatsu InGaAs 
PIN photodiode array (G7150) which actually has 16 photodiodes 
however we choose to use only 8 of these in the present experiment. 
The filling factor for the array is 90%, meaning that 10% of the 
light does not hit an active surface. The quantum efficiency for 
the photodiodes are 80%. 

Virtual Networks. For even numbered mode bases 
(TV = 2,4,6,8) the method for creating the virtual network is as 
follows. The two squeezed modes a\ and 0,2 are combined on a 
half reflecting beamsplitter (HBS). As the output of this HBS is 
an EPR state we choose to call this the EBS. The EBS outputs are 
symmetrically combined with N — 2 vacua, as in Figure 1 of the 

main text. The BSs are then given by ^(cos^ 1 !/ \J ^ — n), where 
n is the number of BSs between the EBS and the BS in question. 
For N=4 and N=6 and N=8, mode output 2 is swapped with 
mode output N-l. For N=8, an additional swap of output modes 
4 and 5 is made. For odd numbered mode bases (N = 3, 5, 7) the 
method is the same with the following modifications. The EBS has 
its reflectivity changed to r = | — . (See for example references 
|29l I30| for more details on N = 3) . The vacua are mixed using 
beamsplitters as above, with one output arm having one less 
vacuum input, n phase shifts are applied to all BS outputs on 
the left of the EBS except for the one left output exiting the last 
BS. Mode outputs 1 and N-l are swapped, and the network for 
N=7 has an additional swap between output modes 3 and 4. The 
homodyne gains g\ are optimised using a genetic algorithm in order 
to maximally satisfy the van Loock-Furusawa inequalities. These 
gains g\ scale the contributions of the quadrature variances and 
are independent from calculations regarding U ne t- Here, optimal 
homodyne gains are calculated using two measures: minimising 
the mean of the N — 1 inequalities; and minimising the variance 
of the set of inequalities. A trade-off between the two measures 
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is needed, and preference is given to minimising the mean of the 
inequalities. 

Spatial mode bases. The input matrix is defined as follows: 



The criteria for verifying the measurements of the various 
cluster states are given below |15II17I with the results summarised 
in Table 2: 
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(5) 



with U* n U- la = I. Each row of U- ln represents the 8 electronic 
gains that match the detection basis to the input modes. For exam- 
ple the top row containing all ones recovers the standard Gaussian 
TEMoo mode (GM), and the second row recovers the phase-flipped 
Gaussian mode (FM). By setting U n et = I we can label each row 
of U- m as G™. Formally, = (In N ^ 8 _ N) )U in , where Nj(8 _ N ) 
is a zero matrix of size N by (8-N). 

The linear optics network for the ideal and symmetric 2-mode 
EPR basis is simply a HBS: 



From equation (2) we get: 



1 1 

V2 s/2 

1 1_ 

V2 \/2 



(6) 
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from which we see that indeed G 2 and G\ indeed share no part 
of the detected light. We show this ideal EPR basis in Figure [3] 
in order to emphasise the spatial separation. Note that the factor 
-^g has been omitted from the scale in Figure j^j for clarity. 

The optimised network uses a beam-splitter reflectivity of 48.8%, 
and produces the following output modes: 



.2. opt 
.2, opt 



— r/ 2 ' opt 77. 2 i 
s/8 nct m 

j_ ( 0.699 0.716 
s/8 I 0.716 -0.699 



11111 1 1 1 
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0.017 
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0.017 



1.414 
0.017 



-0.017 
1.414 



-0.017 
1.414 



The Fourier transform is a rotation of ^ de- 

Therefore we 
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Measuring cluster states. The difference between a two- 
mode cluster state and a two-mode EPR state is a Fourier trans- 
form on one mode. 

grees: F = R {«/2) = {ZZ%- S c2% 2 ) 

get F ( p) = It follows that the homodyne measurements 

we perform in the 2-mode cluster basis {pi — X2,P2 — %l}cluster 
and the 2-mode EPR basis {x± — X2,Pi + P2} EPR are equivalent. 
Therefore we may perform local Fourier transforms so long as we 
can match the homodyne detection basis for individual modes. It is 
important to note that this convenient basis change will not always 
be possible for different clusters. However, by shaping the local os- 
cillator we may have access to arbitrary cluster states within the 
one-beam. This was out of the scope for the current experiment. 



[9 



[10 



-0.017 -0.017 
1.414 1.414 



[11 
[12 

[13 
[14 



[is: 

[16 



N=2 

I([A(pi - £ 2 )] 2 } + <[A(p 2 - £i)] 2 >< 1, 
N=3 

I([A(pj - x 2 )] 2 ) + ([A(p 2 - £1 - x 3 )] 2 ) < 1, 
II([A(p 2 - X! - f 3 )] 2 > + <[A(p 3 - £ 2 )] 2 } < 1, 

N=4 

I([A(pi - £ 2 )] 2 > + GA(p 2 -X! - £ 3 )] 2 > < 1, 
II([A(p2 - xi - £ 3 )] 2 > + ([A(p 3 - £2 - x 4 )] 2 ) < 1, 
III([A(P3 - x 2 - x 4 )} 2 ) + ([A(p 4 - x 3 )] 2 ) < 1, 

N=5 

I([A( P1 - £ 2 )] 2 > + ([A(p 2 - £•! - £ 3 )] 2 > < 1, 
II([A(p2 - £1 - x 3 )] 2 ) + <[A(p 3 - x 2 - x 4 )] 2 ) < 1, 
III([A(p 3 -£2- £ 4 )] 2 > + <[A(p 4 -X3- £ 5 )] 2 > < 1, 
IV([A(p4 -X3- x 5 )} 2 ) + <[A(p 5 - £ 4 )] 2 > < 1. 
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FIG. 1: Multimode entanglement via emulated linear optics networks. Squeezed light and vacua are mixed 
together using unitary operations in order to produce entangled mode states. Unless otherwise stated beam-splitters are 
50% reflective. Superscripts denote mode basis and subscripts denote mode number, (a) The emulated linear optics 
network used to measure 2-mode EPR entanglement (U^ ct ). (b) 8-mode entanglement via a calculated concatenation of 
beam-splitter and tv phase shift operations (f7®et)- The dots between a and b imply virtual networks for N = 3. ..7, not 

shown for brevity. 




FIG. 2: Experimental setup, (not to scale). Squeezed light is prepared and combined in squeezers with a piezo electric 
transducer (PZT) controlling the phase between the two squeezed modes, locked in quadrature. Vacuum modes (vac) 
co-propagate so that the beam exiting squeezers and entering detection contains 8 measurable spatial modes. Multi-pixel 
homodyne detection (MPHD) is used to measure the quadrature amplitudes of the beam in 8 different regions, in detection. 
Local oscillator (LO) gives a reference to phase quadratures. A PC is used to calculate electronic gain functions G„ via the 
notion of virtual networks. The detected beam is then projected onto a basis of measured modes (see equation 1). (Inset) 
Flip mode (FM) generation; half of the wave is phase retarded by half a wavelength, flipping the electric field amplitude. 



Mode pattern matching 




FIG. 3: Spatial mode patterns. Measured modes are defined by spatial patterns of electric field amplitudes. Shown in 
the mode pattern matching box is an example of how the spatial mode pattern for a$ is matched by applying 8 electronic 
gain values (G5) to the detected Gaussian profile (i). The basis of input modes ai...&8 is shown in the middle row (see 
Methods). The arrows represent a mapping via the virtual networks Unct an d U^ct onto the respective bases of entangled 
modes; the top row shows the symmetric EPR or 2-mode basis, while the bottom row shows the 8-mode basis. There is a 
one-to-one correspondence between spatial mode bases shown here ({ai},{af}, {af}) and those shown in Figure]!] Again, 

spatial mode bases for iV=3 to N=7 not shown for brevity. 
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FIG. 4: Noise variance measurements of the spatial modes, (a) x quadrature measurements of the input mode 
basis. The squeezed ([Aa;i] 2 } is shown in the red and anti-squeezed ([A22] 2 } is shown in the blue. The x quadrature 
variances of the 6 vacua modes are measured to equal quantum noise (OdB). (b) p quadrature measurements. The anti- 
squeezed {[Api] 2 } is shown in the red and squeezed ([AP2] 2 } is shown in the blue. The p quadrature variances of the 6 
modes are again measured to equal quantum noise, confirming they are vacua, (c) These variances show the x quadrature 
correlations between modes as in the first half of the L.H.S. of equation (2) of the text. Every column shows N-l traces of x 
quadrature correlations below shot noise, as well as the blue shot noise trace (OdB) normalised to two units of vacua. Each 
green trace shows ([A(a;i — i?)] 2 ) f° r each N-mode basis. Each new colour represents the other N-l variance correlation 
traces of equation (2). (d) Correlations between measured modes in p quadrature, second half of the L.H.S. of equation 



(2). Each green trace now shows ([A(pi +p 2 + gzPz + 



•■+5nPn) 



) . The traces overlapping show that each pair of modes 



is entangled with the same strength as any other pair of modes, a result of optimising for symmetry in the virtual networks. 
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FIG. 5: Inseparability for different entangled mode bases. The solid black line represents the bound of separa- 
bility. Dashed lines represent theory, (a) The blue markers are the averaged measured experimental values for N-mode 
inseparability (right column of Table 1), and the dashed blue line joins the theoretical values of inseparability with the 
same two squeezed inputs used in the experiment. All experimental losses have been taken into account. The red circles 
are the measured experimental values for N-mode cluster states with theory indicated again by the dashed red line. Here 
the maximum value of each row in Table 2 is shown rather than the average value, in order to show that cluster states have 
a much more stringent requirement on squeezing levels (N=5 is clearly separable and not a cluster state), (b) All traces 
have two squeezed inputs and N-2 vacua modes, as in the experiment. What changes is the amount of squeezing in the two 
squeezed inputs, assumed here to be symmetric with equal anti-squeezing. From the top we have: -ldB (magenta); -3dB 
(cyan); experimental parameters (blue); experimental values (blue markers); -6dB (green); and -lOdB (red). 



